/* **************************************************************************
	
	Name of file: 	deal_pp_replication.do
	
	Author: 		Cameron Deal (cameron.h.deal@vanderbilt.edu)
									
	Project: 		Bostock and Public Opinion
	
	Dataset:     	Nationscape 2019-2021
	
	Last Edited: 	24/April/2023
	
	Note:			
	
*************************************************************************** */

capture log close
clear all
set more off
clear mata
clear matrix
ssc install estout, replace
ssc install coefplot, replace
set matsize 11000

*set directory for data
cd "/*INSERT YOUR DIRECTORY HERE*/"

*load data
use "Data\nationscape.dta", clear

*set results directory 
cd "/*INSERT YOUR DIRECTORY HERE*/"

*set style graphs
set scheme /*sj*/ /*s1mono*/ s2mono

************************************************
* Data Cleaning
************************************************

*****Outcomes

*Alternative Outcome Variable (continuous)
gen lgbt_fav_cont=.
replace lgbt_fav_cont=0 if group_favorability_lgbt==4 /*very unfavorable*/
replace lgbt_fav_cont=1 if group_favorability_lgbt==3 /*somewhat unfavorable*/
replace lgbt_fav_cont=2 if group_favorability_lgbt==2 /*somewhat favorable*/
replace lgbt_fav_cont=3 if group_favorability_lgbt==1 /*very favorable*/
label var lgbt_fav_cont "LGBT Favorability Continuous"

*Main Outcome Variable- Stigma towards LGBT People (Binary) omitting those who don't answer the question or say "I haven't heard enough"
gen lgbt_fav_bin=.
replace lgbt_fav_bin=0 if lgbt_fav_cont==1 | lgbt_fav_cont==2 | lgbt_fav_cont==3 /*all else*/ 
replace lgbt_fav_bin=1 if lgbt_fav_cont==0  /*very unfavorable*/
label var lgbt_fav_bin "LGBT Stigma"

*****Treatment/Post Variables

*Post Bostock (Binary)- Right at midnight on June 14th
egen response_id_num = group(response_id)
gen post=0
replace post=1 if response_id_num>=291235 

*Treatment (Binary)- States treated by Bostock are 1, those who already had protections are 0.
gen treatment= 0
replace treatment = 1 if inlist(state, "AK","MT","ND","MI","ID","WY","SD")
replace treatment = 1 if inlist(state,"IN","OH","PA","NE","MO","KY","WV")
replace treatment = 1 if inlist(state,"KS","AR","TN","NC","SC","AZ")
replace treatment = 1 if inlist(state,"OK","LA","MS","AL","GA","TX","FL")

*Difference-in-Differences Estimator
gen did= post*treatment

*State- Recoded to factor variable
encode state, gen(state_factor)

*Week- Formatted in week_formatted
tostring week, generate(week_1)
replace week_1="20190801" if week_1=="20190800"
replace week_1="20201001" if week_1=="20201000" 
replace week_1="20191031"  if week_1=="20191032"  
gen week_formatted= date(week_1, "YMD")
format week_formatted %td

*Trend Variable (Starts at 1 in Week 1 and continues to 77 in Week 77)
gen int trend = wofd(week_formatted) - wofd(td(11jul2019))
replace trend = 77 if trend==78

*****Individual Controls

*Age- Continuous is already in the dataset, generate age squared.
gen age_squared= age^2

*Unemployment (Binary)
gen unemployed= (employment==4)

*Education- 3 Categories (High School or Less, Some College, and Bachelor's Degree or Higher)
gen educat=.
replace educat = 1 if education <= 4
replace educat = 2 if education <= 7 & education>=5
replace educat = 3 if education>=8
label define educat 1 "High School or Less" 2 "Some College" 3 "Bachelor's or Higher"
label values educat educat
label var educat "Education Categories"

*Race/Ethnicity- 4 Categories (White NH, Black NH, Other Races NH, and Hispanic)
gen race4 =.
replace race4 = 1 if race_ethnicity==1 & hispanic==1
replace race4 = 2 if race_ethnicity==2 & hispanic==1
replace race4 = 3 if race_ethnicity>=3 & hispanic==1
replace race4 = 4 if hispanic>=2
label define race4 1 "White, Non-Hispanic" 2 "Black, Non-Hispanic" 3 "Other Races, Non Hispanic" 4 "Hispanic"
label values race4 race4
label var race4 "Racial Categories"

*Religion- 5 categories (Protestant, Catholic, Other Christian, Other Religions, No Religion)
gen relig5 =.
replace relig5 = 1 if religion==1
replace relig5 = 2 if religion==2
replace relig5 = 3 if religion>=3 & religion<=5
replace relig5 = 4 if religion>=6 & religion<=9
replace relig5 = 4 if religion==13
replace relig5 = 5 if religion>=10 & religion<=12
label define relig5 1 "Protestant" 2 "Catholic" 3 "Other Christian" 4 "Other Religions" 5 "No Religion"
label values relig5 relig5
label var relig5 "Religious Categories"

*Interest in government (Binary)
gen interested = (interest==1) /*they follow what's going on with the government most of the time*/

*Party Identification (Binary Democrat and Republican- Reference group is other political party+independents+something else)
gen democrat = (pid3==1)
gen republican = (pid3==2)

*Combine weeks
egen weeks_bin3 = cut(trend), at(0,3,6,9,12,15,18,21,24,27,30,33,36,39,42,45,48,51,54,57,60,63,66,69,72,75,78) icodes


************************************************
* Global
************************************************

*Individual-Level Controls
global indiv_controls "i.gender age age_squared i.unemployed i.educat i.race4 i.relig5 i.democrat i.republican"

************************************************
* Restrict to complete cases and variables of interest
************************************************
*Variables of Interest 
keep lgbt_fav_bin lgbt_fav_cont did trend week_formatted week state_factor age age_squared unemployed gender educat race4 hispanic relig5 weight census_region post treatment interested democrat republican weeks_bin3 state congress_district
			
			
*Complete Cases
regress lgbt_fav_bin i.post i.treatment did $indiv_controls i.state_factor i.week [pweight=weight], vce(cluster state)
	gen s = 1 if e(sample)
	drop if s==.
	
**********************************************************************************************************************************************************************************
**********************************************************************************************************************************************************************************
**********************************************************************************************************************************************************************************
**********************************************************************************************************************************************************************************
**********************************************************************************************************************************************************************************
************************************************
* Figures 1 and 2: Distribution of Disapproval
************************************************
*Figure 1: Treatment vs. Untreated
preserve
collapse (mean) lgbt_fav_bin treatment post, by(state_factor)
kdensity lgbt_fav_bin, nograph generate (x fx)
kdensity lgbt_fav_bin if treatment==0, nograph generate(fx0) at(x)
kdensity lgbt_fav_bin if treatment==1, nograph generate(fx1) at(x)
label var fx0 "Untreated"
label var fx1 "Treated"
line fx0 fx1 x, sort ytitle(Density) xtitle(% Very Unfavorable) graphregion(color(white)) bgcolor(white) legend(region(lwidth(none))) 
graph export fig1.png, replace
restore

*Figure 2a: Before vs. After- Treated
preserve
collapse (mean) lgbt_fav_bin treatment, by(state_factor post)
kdensity lgbt_fav_bin, nograph generate (x fx)
kdensity lgbt_fav_bin if treatment==1 & post==0, nograph generate(fx0) at(x)
kdensity lgbt_fav_bin if treatment==1 & post==1, nograph generate(fx1) at(x)
label var fx0 "Before"
label var fx1 "After"
line fx0 fx1 x, sort ytitle(Density)  xtitle(% Very Unfavorable) graphregion(color(white)) bgcolor(white) legend(region(lwidth(none))) 
graph export treated.png, replace
restore

*Figure 2b: Before vs. After- Untreated
preserve
collapse (mean) lgbt_fav_bin treatment, by(state_factor post)
kdensity lgbt_fav_bin, nograph generate (x fx)
kdensity lgbt_fav_bin if treatment==0 & post==0, nograph generate(fx0) at(x)
kdensity lgbt_fav_bin if treatment==0 & post==1, nograph generate(fx1) at(x)
label var fx0 "Before"
label var fx1 "After"
line fx0 fx1 x, sort ytitle(Density)  xtitle(% Very Unfavorable) graphregion(color(white)) bgcolor(white) legend(region(lwidth(none))) 
graph export untreated.png, replace
restore

**********************************************************************************************************************************************************************************
**********************************************************************************************************************************************************************************
**********************************************************************************************************************************************************************************
**********************************************************************************************************************************************************************************
**********************************************************************************************************************************************************************************
************************************************
* Figure 3: Heterogeneous Effects Across Observables
************************************************
*Government Interest
regress lgbt_fav_bin i.post i.treatment did $indiv_controls i.state_factor##c.trend [pweight=weight] if interested==1, vce(cluster state)
eststo gov_interest
regress lgbt_fav_bin i.post i.treatment did $indiv_controls i.state_factor##c.trend [pweight=weight] if interested==0, vce(cluster state)
eststo notgov_interest
*esttab using table2reg_022222.csv, se keep(did) r2 star(* 0.10 ** 0.05 *** 0.01) label                 ///
 *   title(Government Interest)  
*eststo clear
mean lgbt_fav_bin if post==0 [pweight=weight], over(interested)

*Party Affiliation
regress lgbt_fav_bin i.post i.treatment did $indiv_controls i.state_factor##c.trend [pweight=weight] if democrat==1, vce(cluster state)
eststo democrat
regress lgbt_fav_bin i.post i.treatment did $indiv_controls i.state_factor##c.trend [pweight=weight] if republican==1, vce(cluster state)
eststo republican
*esttab using table2reg_022222.csv, se keep(did) r2 star(* 0.10 ** 0.05 *** 0.01) label                 ///
 *   title(Party Affiliation)   append
	*eststo clear
mean lgbt_fav_bin if post==0 & democrat==1 [pweight=weight]
mean lgbt_fav_bin if post==0 & republican==1 [pweight=weight]

*Gender 
regress lgbt_fav_bin i.post i.treatment did $indiv_controls i.state_factor##c.trend [pweight=weight] if gender==1, vce(cluster state)
eststo female
regress lgbt_fav_bin i.post i.treatment did $indiv_controls i.state_factor##c.trend [pweight=weight] if gender==2, vce(cluster state)
eststo male
*esttab using table2reg_022222.csv, se keep(did) r2 star(* 0.10 ** 0.05 *** 0.01) label                 ///
 *   title(Gender)   append
	*eststo clear
mean lgbt_fav_bin if post==0 [pweight=weight], over(gender)

*Race=White
regress lgbt_fav_bin i.post i.treatment did $indiv_controls i.state_factor##c.trend [pweight=weight] if race4==1, vce(cluster state)
eststo white
regress lgbt_fav_bin i.post i.treatment did $indiv_controls i.state_factor##c.trend [pweight=weight] if race4!=1, vce(cluster state)
eststo nonwhite
*esttab using table2reg_022222.csv, se keep(did) r2 star(* 0.10 ** 0.05 *** 0.01) label                 ///
 *   title(Race) append
	*eststo clear
mean lgbt_fav_bin if post==0 & race4==1 [pweight=weight]
mean lgbt_fav_bin if post==0 & race4!=1 [pweight=weight]

*Education
regress lgbt_fav_bin i.post i.treatment did $indiv_controls i.state_factor##c.trend [pweight=weight] if educat==3, vce(cluster state)
eststo bachelors
regress lgbt_fav_bin i.post i.treatment did $indiv_controls i.state_factor##c.trend [pweight=weight] if educat!=3, vce(cluster state)
eststo no_bachelors
*esttab using table2reg_022222.csv, se keep(did) r2 star(* 0.10 ** 0.05 *** 0.01) label                 ///
*    title(Education) append
*	eststo clear
mean lgbt_fav_bin if post==0 & educat==3 [pweight=weight]
mean lgbt_fav_bin if post==0 & educat!=3 [pweight=weight]

*Above/Below Median 
gen below_median=.
replace below_median = 0 if treatment==0 
replace below_median = 1 if inlist(state_factor, 1, 4, 10, 11, 14, 16, 17, 23, 29, 30, 36, 39, 44, 51)
gen above_median=.
replace above_median = 0 if treatment==0 
replace above_median = 1 if inlist(state_factor, 2, 3, 18, 19, 25, 26, 27, 28, 37, 41, 42, 43, 50)

regress lgbt_fav_bin i.post i.treatment did $indiv_controls i.state_factor##c.trend [pweight=weight] if below_median!=., vce(cluster state)
eststo below_median
regress lgbt_fav_bin i.post i.treatment did $indiv_controls i.state_factor##c.trend [pweight=weight] if above_median!=., vce(cluster state)
eststo above_median

label var did "DiD"

*Coefficient Plot
coefplot (gov_interest, label(Interested) pstyle(p1))  ///
(notgov_interest, label(Not Interested) pstyle(p1))  ///
(democrat, label(Democrat) pstyle(p2))  ///
(republican, label(Republican) pstyle(p2))  ///
(female, label(Female) pstyle(p3))  ///
(male, label(Male) pstyle(p3))  ///
(white, label(White) pstyle(p4))  ///
(nonwhite, label(Non-White) pstyle(p4))  ///
(bachelors, label(Bachelors) pstyle(p5))  ///
(no_bachelors, label(<Bachelors) pstyle(p5))  ///
(above_median, label(Above Median) pstyle(p6))  ///
(below_median, label(Below Median) pstyle(p6))  ///
, keep(did) xline(0) graphregion(color(white)) bgcolor(white) legend(region(lwidth(none)) position(9) cols(1) rowgap(3.3pt)) 
graph export heterogeneity.png, replace
******************************************************************************************************************